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Abstract 


A knowledge of thermodynamic properties is of immense importance for the 
design and development of processes and equipment. In several situations in which no 
experimental data for the necessary properties is at hand, the values must be predicted. 
The work reported in this thesis is an attempt to develop a computer software to 
estimate thermodynamic properties of importance in chemical engineering, using cubic 
equation of state. The cubic equations of state used for the estimation of properties are 
van der Waals equation, Redlich Kwong equation, Soave Redlich Kwong equation and 
Peng Robinson equation. The software includes programs for predicting properties for 
pure components as well as multi-component mixtures. For pure substances the 
properties included are molar volume, departure functions, fugacity, saturation 
pressure, saturation temperature, vapor enthalpy and vapor entropy. For multi- 
component mixtures, it includes molar volume, departure function, fugacity and partial 
fugacity. In addition to performing the phase equilibrium calculations the software is 
also capable of predicting vapor-liquid equilibrium data (VLE). For the VLB data 
prediction it uses five different activity coefficient models including UNIFAC. The 
software also contains an inbuilt database, which includes 250 compounds. 
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Chapter 1 

INTRODUCTION 


The main activities of a chemical engineer are concerned with the design and 
development of new processes or the improvement of existing processes for converting 
the raw materials into useful products of greater value in an optimal way. The resources 
available to an engineer are the raw materials, space and time. The emphasis should be 
given on consuming minimum amount of resources for producing the desired products. 
But before initiating the task one should know its viability, the energy requirement and 
the fractional conversion of the raw materials into product, for the process. For example, 
designing a process for producing alcohol by a specific reaction requires the knowledge 
of the feasibility of the reaction at a particular temperature and pressure, the extent to 
which the reaction can proceed, the energy requirement for the reaction, the effect of 
change in operating condition and energy required to separate the desired product from 
the byproducts and the reactants. All these information can be obtained by the 
application of the principles of thermodynamics. 

In thermodynamic analysis of a process, it is essential to have knowledge of the 
properties of the system and the relationships that exist among the various properties. 
The properties of interest in thermodynamics can be broadly classified into two groups 
as 

> Measurable Properties and 

> Non-measurable Properties. 

The properties which can be directly measured, like pressure (P), temperature 
(7) and Volume (V), are called measurable properties and the properties, that cannot be 
directly measured, like internal energy (U), enthalpy (B), Helmholtz free energy (A), 
entropy (S), Gibbs free energy (G), etc., are called non-measurable properties. The 
properties are independent of the path followed by a system in reaching a given state as 
the change in the properties depends only on the initial and the final state of the system. 
As a property has a unique value when a system is in a given state, the non-measurable 
properties can be estimated in terms of the measurable properties, T, v and P. There 



exist fundamental property relations, which contain all the information about the system 
and if they are known, in principle, each and every property of the system can be 
estimated. There are several alternative fundamental relations. The most commonly 
used fundamental relations are 


U=U(S,V,Ni,N2.Ns, . 

... NrJ 

(1.1) 

H = H (S. P, Ni, N 2 , Ns, . 

.■,Nr) 

(1.2) 

A =A(T,V,Ni,N2.N3, .. 


(1.3) 

G-G(T, P, Nu N2, Ns. . 

....Nr) 

(1.4) 

S = S(U,V.NuN2.Ns, .. 

.,N^ 

(1.5) 


Where, 

N„ is the number of moles of the component n in the system. 

The first relation is the fundamental relation in the energy representation and the next 
three are the partial Legendre transforms of the first relation. The last relation is the 
fundamental relation in the entropy representation. First four relations are called 
thermodynamic potentials and they can be effectively used for developing the 
relationships among the thermodynamic properties. For a simple compressible system 
H, S and V can be written as a function of T and P as 



H = H(T,P). 

(1.6) 


S (T. P) 

(1.7) 

and 

V= V(T,P) 

(1.8) 


If H, S and V are known for a particular system, then, G, U and A can be 
calculated from the following relations 


1 

II 

(1.9) 

U=H-PV 

(1.10) 

A = U-TS 

(1.11) 


Usually, the thermodynamic property data is presented in the following three forms 


> ' Equation Of State (EOS), 

> Diagrams and charts and 

> Thermodynamic property tables. 



The EOS axe the relationships between P, v and T, having the form/ (P, v, T) 
=0, and are valid over a wide range of temperature and pressure. Several forms of EOS 
are available and they are discussed in Chciptsf 2. There are several thermodynamic 
charts, tables and diagrams that provide the values of the properties for various ranges 
of temperature and pressure. For example the Mollier diagram H vs S curves using 

T and P as a parameters; steam tables provide the properties like saturation pressure, 
specific volumes, specific enthalpy and specific entropy for the steam. The estimation of 
thermodynamic properties using the data available in the form of charts, tables and 
diagrams requires frequent interpolation. Moreover, the estimation of some properties, 
which are not presented in the charts and tables, requires knowledge of the derivatives 
of some other properties, which cannot be evaluated with reasonable accuracy. On the 
other hand EOS can be directly used for the calculation of various thermodynamic 
properties. Moreover some EOS are capable of expressing liquid as well as the vapor 
phase behaviors and use only three compound dependent parameters like critical 
temperature (Tc), critical pressure {Pc) and acentric factor {d). 

In the present work, attempts have been made to develop a computer software 
for the prediction of the thermodynamic properties using four pressure explicit cubic 
EOS. The software is developed using Borland C++ programming language on MS- 
DOS operating system. It can predict the thermodynamic properties for pure substances 
as well as for multi-component mixtures. For the pure components it can predict molar 
volume, departure functions, fiigacity, saturation pressure, saturation temperature, vapor 
enthalpy and vapor entropy. For multi-component mixtures it can predict molar volume, 
departure functions, fugacity, fiigacity of the individual components in the mixture 
(partial fugacity). In addition to performing phase equilibrium calculations that is 
bubble point, and dew point calculation the software also predicts the constant pressure 
Vapor-Liquid Equilibrium (VLE) data for binary systems using various activity 
coefficient models. The software provides an option to choose anyone of the following 
four cubic EOS. 

> van der Waals equation of state 

> Rediich-Kwong equation of state 

> Soave- Redlich-Kwong equation of state and 



> Peng-Robinson equation of state 


For the predication of constant pressure VLB data the software provides the following 
options for the activity coefficient models. 

> Margules equation 

> van Laar equation 

> Wilson equation and 

^ Non-Random Two Liquid (NRTL) equation and 

> UNIFAC model (A group contribution model) 

The software has an inbuilt database for 250 compounds and can process the relevant 
information to display the desired thermodynamic property at specified temperature and 
pressure. If one selects a compound that is not included in the database, the user has to 
furnish the compound related parameters. For the constant pressure VLB data using 
activity coefficient models the user have to provide the parameter related to the model 
chosen. However, for the UNIFAC model, a separate database containing more then 70 
compounds is included in the software. 



Chapter 2 

THEORY 


For the estimation of thermodynamic properties various theoretical and 
empirical relations, although not generally valid, are in common use. Most of the better 
estimation methods use equations based on theory coupled with empirical correlation of 
the constants that are not provided by the theory. Introduction of empiricism into parts 
of a theoretical relation provides a powerful method for developing a reliable 
correlation. The best example is the various equations of state, which are most widely 
used for the prediction of thermodynamic properties. The main advantages of using an 
equation of state are already discussed in the previous chapter. 

The first section of this chapter deals with the various equations of state 
followed by the estimation methods for thermodynamic properties. 

2.1 Equation Of State: 

The most convenient method of representing the P-v-T behavior of a substance 
is through a mathematical expression of the type 

f(P, v,T) = 0 (2.1) 

This type of relation is called Equation Of State (EOS). The simplest and the oldest 
volumetric EOS is the ideal gas equation, given by 

PV = nRT (2-2) 

where, R is the universal gas constant. 

The Eq. (2.2) is based on the assumptions that volume occupied by the 
molecules themselves is negligible, molecules are independent of each other and do not 
exert any attractive force. The ideal gas equation is valid only at sufficiently low 



density. Moreover, it does not predict the occurrence of a vapor liquid phase transition. 
As the ideal gas law is only an approximation or idealization to the behavior of the 
gases and is applicable only at low densities, several attempts have been made to 
accurately represent the behavior of real gases and a large number of EOS have been 
proposed. Following are the most commonly used cubic EOS. 

> Van der Waals Equation of state (vdW) 

> Redlich-Kwong Equation of state (RK) 

> Soave- Redlich-Kwong Equation of state (SRK) 

> Peng-Robinson Equation of state (PR) 

2.1.1 van der Waals Equation of state 


The van der Waals EOS, proposed in 1873, was the first EOS capable of 
representing the vapor liquid coexistence. It is given by 


P=IL.± 

v-b 


(2.3) 


The EOS parameter a is a measure of attractive forces between the molecules, 
and the parameter b is the excluded volume i.e. the volume that is not available for the 
molecular motion due to the presence of the other molecule. The parameters a and can 
be detemined using the fact that the critical isotherm exhibits a point of inflection. That 
is, at the critical point 


\dv 


= 0 and 




Jr, 


= 0 


(2.4) 


Using Eq. (2.3) and Eq. (2.4) one can obtain 


IIPPT} 
64Pr 


a 


(2.5) 


and 




iPr 


( 2 . 6 ) 



where, 7c and Pc are the critical temperature and pressure, respectively. 


2.1.2 Redlich-Kwong Equation of state 

An important modification of the vdW EOS was made by Redlich and Kwong 
(1949). The RK EOS retains the original repulsive term with the addition of a 
temperature dependent attractive term. The RK EOS is given as 


L_ 

v-b v{v + b) 

Where 


a = 0.42748 ^ 

Pc 

JDT 

b = 0.08664 — ^ 


(2.7) 


( 2 . 8 ) 

(2.9) 


2.1.3 Soave- Redlich-Kwong Equation of state 


Soave (1972) proposed a modification to the RK EOS replacing the term — 
with a more general term a (T), that is 

P = ( 2 . 10 ) 

v-b v(v + ^>) 

Where 

a{T) = a,a(T) (2.11) 

a(r) = [l + (o.48 + 1.574®-0.176(a^)(l-V^)] (2-12) 

And ty is the Pitzer (1977) acentric factor, which is defined as 

£n = — lOgP,. I 7,= 0.7 ^ 


(2.13) 



where, P/ (reduced vapor pressure at 7; = 0. 7) = 


The SRK EOS predicts the vapor pressure of a large number of hydrocarbons, 
which are fairly close to the experimental data. 

2.1.4 Peng-Robinson Equation Of State 


Peng and Robinson (1976) used a different volume dependence to give slightly 
improved liquid volume (Zc=0.30T) and redefined the term a (T). The PR EOS is given 
by 


Where, 


p_ RT a{T) 

v-b v(v + Z>)+Z>(v-6) 

(2.14) 

R%^a 

(2.15) 

a = 0.45724 ^ — 

Pc 

RT 

b = 0.07780 ^ 

(2.16) 

Pq 


(2.17) 

S = 0.37464 + 1 .5426iy - 0.269fi)' 

(2.18) 


2.2 Calculation of thermodynamic properties 
for pure substances 

The four cubic EOS, discussed in the previous section can be expressed in a 
general form as 

p ^ AL - - (2.19) 

(v-b) (v^ + ubv + wb^) 

Where, m and w are constants associated with a given EOS. The values of u and 
w for the four EOS are presented in Table 2.1. 



The Eq. (2. 19) can be rewritten as 


Where, 


+aZ'^-\- p Z ■¥ y= 0 


Z 

a 

P 

7 

A 

B 


RT 

= -I -B + uB 
= A + coB^ -uB-uB^ 
= -AB-coB'^-coB^ 

aP 

(RTf 

bP 

RT 


( 2 . 20 ) 

( 2 . 21 ) 

( 2 . 22 ) 

(2.23) 

(2.24) 

(2.25) 

(2.26) 


Table 1.0 Values ofu,w, cx, fi and y cubic equations of state 


Equntion of stnte 

u 

w 


P 

y 

Van der Waals 

0 

0 


A 

-AB 

Redlich-Kwong 

1 

0 

-1 


-AB 

Soave-Redlich-kwong 

1 

0 

-1 

EBUIH 

-AB 

Peng-Robinson 

2 

-1 


gglglll 



The Eq. (2.20) is the general form of the four cubic EOS. The calculation 
method presented in the following sub-section uses this general form of EOS. 

2.2.1 Calculation of molar volume 

The volume of the liquid or vapor can be calculated by adopting an iterative 
technique. However, such calculations are tedious and hence one can use Cardan s 
method to determine the roots of cubic EOS. 

The Eq. (2.20) can be solved using Cardan’s method and depending upon the 

pressure and temperature one can get 






















Only one real and positive root. If 7^ Tc it corresponds to vapor 
compressibility factor Zv else it will corresponds to either liquid or vapor 
compressibility factor. 

■Z Three real and positive roots, two of them are equal. The smallest value 
of Z corresponds to Z/ and the largest value corresponds to Zy. 

■Z Three real and unequal roots. The smallest value of Z corresponds to Z; 
and the largest value corresponds to Zy. The intermediate value has no 
physical significance. 

These values of Zy and Z; are then used to predict vapor and liquid volumes 
respectively using Eq. (2.2 1 ). 

2.2.2 Calculation of Enthalpy and Entropy Departure 
Functions 


The enthalpy departure, which is the difference between the enthalpy of a real 
fluid and the enthalpy of the ideal fluid at the same temperature and pressure, is given 
by 


=i?7’(Z-l) + 


V 



\feoo 



dv 


T 


{ 121 ) 


Similarly, the entropy departure is the difference between the entropy of the real 
fluid and the entropy of the ideal fluid at the same temperature and pressure and for the 
pressure explicit cubic EOS it is given by 




V=oo 



(2.28) 


Following are the expressions for the enthalpy departure and entropy departure 
functions using the four cubic EOS. 



vdW EOS 


RKEOS 


SRK EOS 


where 


PR EOS 


where, 


(h-h\p=Pv-RT-- 

V 

RT 


2b 


(s-s\p=Rln(Z-B)- 


a 


2bT 


•In 


2 ; 

^ Z+B ^ 

Z 


T 


{h-h\p =RT{Z-\) + - 




-a 


-In 


(¥' 


=i?ln(Z-5)-f 


\( da^ 

~b[^) 


In 


(Z + B\ 


\ ^ J 


fda^ 

0A2748R^T^ 1 

[dTj 

Pc i 


TTr 


(2.29) 

(2.30) 

(2.31) 

(2.32) 

(2.33) 

(2.34) 

(2.35) 




T 


fda\ 

^dT) 

24Zb 



V 


z+m+ 4 ^ 

z+bIi-S 


h 

ij 


=i?ln(Z-5)+ 


1 r da ^ 

2j^\dT^ 


In 


^z +5(1-72), 


(2.36) 


(2.37) 


fda^ 

1 0.457245^ Tc [ 

[dT) 

1 " p, i 


a 

TTc 


(2.38) 


The enthalpy and entropy departures can be calculated using the above 
equations. 



2.2.3 Calculation of Fugacity and Fugacity Coefficient 


For an ideal gas, the change in Gibbs free energy is given by 
RT 

dg'" =—dP = RTd\nP 


(2.39) 


Eq. (2.39) is not valid for a real gas. However, its form is quite useful and thus a 
new parameter fiigacity is defined such that the form of the Eq.(2.39) is retained for real 
fluids. Therefore, the fugacity is defined as 

dg = RTd\xif (2.40) 

f 

and p ^ &sP -i^O (2.41) 


The ratio of the fugacity to the pressure is defined as the fugacity coefficient ^ 

that is 

(j)=^ (2.42) 

P 


To estimate the fugacity, from the above discussed four cubic EOS, the 
following relation between the fugacity coefficient, Gibbs free energy and the departure 
function can be used. 


ln(j) = 


g-g 

RT 


h-h’^ Tjs-s’^) 
RT RT 


(2.43) 


By substituting the enthalpy and entropy departure presented in the previous 
section, Eq. (2.29) to Eq. (2.38), in the Eq. (2.43) the following relations for the fugacity 
coefficient can be obtained. 
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vdW EOS 


RKEOS 


SRK EOS 


PR EOS 


RT RTv 


ln,p = Z-l-\n(Z-B ) — ^In 
^ ’ hRT 


^Z+B^ 


\ ^ J 


\n<l> = Z-\-]n(Z-B)-—\n[^^ 

bRT V Z y 


ln{z5 = Z-l-ln(Z-5) — In 


^Z + 5(l + V2^^ 


2^f2bRT [z + b(i-^) 


(2.44) 


(2.45) 


(2.46) 


(2.47) 


2.2.4 Calculation of Saturation Temperature and 
Saturation Pressure 


The criterion for equilibrium can be used for the calculation of saturation pressure and 
saturation temperature of a liquid. If a liquid is in equilibrium with its own vapor, the 
criterion of equilibrium gives 

f = r (2.48) 

Where, 

/^and /'"are the fugacity of the liquid and vapor phases respectively. 

An iterative technique can be used to calculate the saturation pressure (or 
saturation temperature). For the calculation of saturation pressure an initial guess value 
of P is assumed at the specified temperature T (< Tc). Then at these P and T values, 
depending upon the choice of EOS, and /''can be calculated using the Eq. (2.44) to 
Eq. (2.47) discussed in the previous section. Subsequently, it is checked whether the 
criterion of equilibrium specified by Eq. (2.48) is satisfied or not, if the criterion is 
satisfied then the guess value P is the required saturation pressure. If it is not satisfied, 
the guess value will be modified using following equation. 


P{new) = P 


r 


(2.49) 



Similar procedure can be employed to calculate saturation temperature. 


2.2.5 Calculation of Enthalpy and Entropy of Vapor 


For the calculation of enthalpy and entropy the four step path shown in the figure below 
can be used. 



Fig. 2.1 Four-step path for the estimation ofh and s. 


The reference state is assumed to be saturated liquid at 298.15 Kelvin (25 ^’C). At 
this reference state the internal energy and entropy of the substance are assumed to be 
zero. Then 

h‘’=u°+P^^‘v (2.50) 

Where, 

u° = Internal energy at reference state. 

= Saturation pressure at 298.15 K 
V = Molar volume at 298. 15 K and 

which gives 

( 2 - 51 ) 

In the step 1 the substance is vaporized at 298. 15 K and at saturation pressure so 
change in enthalpy is equal to the enthalpy of vaporization at 298.15 K.To calculate the 
enthalpy of vaporization (Ah^) at 298.15 K initially, the enthalpy of vaporization is 
calculated, at nomal boiling point of the substance, using Riedel correlation given by 








AA,„ = 1.093i?r 


(2.52) 


0.930-7;^ 


Then the Watson’s correlation given below is used to determine Ahv. 


AA, = AA., = AA., 


r^_T 

^ ^r2 

vl-^wy 


(2.53) 


The change in entropy can then be calculated by the following relation 


1 j. 


(2,54) 


In the step 2 the substance is brought to the ideal gas state from the real gas state 
at the same temperature and pressure. The changes in enthalpy and entropy can be 
estimated from the knowledge of departure functions. 

In the step 3 the substance is brought to the desired T and P in ideal gas state and 
the changes can be estimated using the constant pressure heat capacity data. The 
following equations can be used for the calculation of the entropy and enthalpy changes 
in the step 3. 


Ah, 


T2 


CldT 


Tl 

T2 


=1 


= \CldT-Rln 


n 


Pr 


(2.55) 

(2.56) 


In the step 4 the substance in the ideal gas state at T and P is changed to the real 
gas sate at the same temperature and pressure. And the changes in the enthalpy and 
entropy can be calculated from the knowledge of departure functions. 

Finally the enthalpy and entropy are obtained by summing all the enthalpy and 
entropy changes in the step 1 to step 4. 


h - ho+ hi + h 2 + hs +h4 
S = Si+S2'^S3+S4 


(2.57) 

(2.58) 



2.3 Calculation of Thermodynamic 
Properties of Mixture 

The thermodynamic properties of the mixture can be obtained by modifying the 
method used for a real gas The modification needs additional composition 

variable as an input parameter. The constants that characterize the pure fluid are 
replaced by the mixture constants. The relation which expresses the mixture constants in 
terms of pure fluid constants are called mixing rules. Usually the mixing mles are 
empirical and their validity rests upon the agreement between predicted and 
experimentally measured properties. The most commonly used mixing rule for the cubic 
EOS (vdW, RK, SRK and PR) is the van der Waals mixing rule. 

In the van der Waals mixing mle the mixture parameters a„ and bm are 
expressed as 

*. =EM ( 2 . 60 ) 

I 

where Ki./is the binary interaction parameter between the components i and j. 

The binary interaction parameter Kij is introduced to obtain better agreement 
between the predicted and experimental results. Usually this parameter is obtained by 
fitting the equation of state to the experimental data. The interaction parameter Kij 
assumes values close to zero and in the absence of any data it can be set to zero. 

The mixture molar volume, departure functions and fugacity can be calculated 
using the methods disused in section 2.2. The prediction of the partial fugacity is 
discussed in the following sub-section. 

2.3.1 Calculation of Partial Fugacity 


The partial fugacity that is the fugacity of component / in the mixture is 
defined as 



~ '^idP - RTd In^. at constant T 


(2.61) 


-^,°)=-^rin A =RT\J-^, 




(2.62) 


1 as P -> 0 


(2.63) 


where, is the fugacity of the pure component i in the reference state, where it 
behaves like an ideal gas. 

The partial fugacity of component i in a mixture for the four cubic EOS can be 
estimated by the relations given below. 

vdW EOS 


RKEOS 


, r b, , iJaUi 

In,^ =— i In Z 1 — — ^ ^ 

v-b I vi RTv 


W,.4(Z-l)-ln(Z-.K^|-2^1nr^ 


(2.64) 


(2.65) 


SRK EOS 


lni=h.(z-i)-ln(Z-B)+— Inf— 

' ^ ^ ’ bRT b \a Z 


PR EOS 


' 4 2-JlRT b 1/a Z + Bl-s/2 


( 2 . 66 ) 


(2.67) 


where 


b. KPa) 

^ ^ yjPj 

i V ^Ci 


( 2 . 68 ) 


To estimate the partial fugacity of a component in a mixture, first the pure component 
constants a,- and bi are calculated. Then using the vdW-mixing rule the mixture 


17 



parameters are calculated. The volume, compressibility factor and h are calculated 
using Eq. (2.64) to Eq. (2.68) depending upon the choice of EOS. 


2.4 Phase Equilibrium Calculation 


For a closed system consisting of co-existing vapor and liquid phases, each 
phase containing C components, the criterion for equilibrium between the two phases is 
given by 

[/ = 1,2,3, ,C] (2.69) 

Using the concept of fogacity, fugacity coefficient, activity and activity 
coefficient, Eq. (2.69) can be transformed to give 

[i = 1,2,3, ,C] (2.70) 

where, 

= Fugacity coefficient of component / in the vapor phase solution and 
= Fugacity coefficient of component i at saturation pressure. 

Eq. (2.69) is called the basic equation of vapor liquid equilibrium and provides a 
systematic means of analyzing vapor-liquid equilibrium data. The analysis of vapor- 
liquid equilibrium data can be divided into the following four problems. 


Problem Type 

Known 

To estimate 

> Bubble pressure calculation 

T, x, ’s 

P,yis 

> Bubble temperature calculation 

P, XiS 

T,y/s 

> Dew pressure calculation 

T.yiS 

P, X/’s 

> Dew temperature calculation 

P,yYs 

T, x,’s 



2.4.1 Bubble point and Dew point Calculation using 
EOS 


For bubble point and dew point calculation basic equation of vapor liquid 
equilibrium, Eq. (2.70), can be rewritten as 


(f/x,. = ^/y. [i = 1,2,3, ,C] 

Or K,=^ = i- 

X, 

In addition to this, the values of x,- andj; must satisfy the stoichiometric conditions 
given by 


c 



1=1 

(2.73) 


c 


and 

1 = 1 

(2.74) 


(j)/ and can be calculated using equations discussed in section 2. 3. 7. An 
iterative scheme can be used using the Eq. (2.72) for the bubble point and dew point 
calculation. 


(2.71) 

(2.72) 


2.4.2 Prediction of VLE data - Excess Gibbs Free 
Energy Model 

The VLE data can be predicted by using Eq. (2.70). At low to moderate pressure 
the Poynting correction factor can be assumed as unity and (<^//<^'/“!) is approximately 
equal to unity. Hence, Eq. (2.70) takes the following form 
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The Eq. (2.75) can be used for the prediction of constant pressure VLE data using 
various excess Gibbs free energy models discussed in the following sub section. 


2.4.2. 1 Excess Gibbs Free Energy Models 

The behavior of almost all liquid solutions cannot be adequately expressed in 
terms of an equation of state. For such mixtures, attempts have been made to estimate 
the excess Gibbs free energy either empirically or semi empirically leading to several 
models for excess Gibbs free energy. In this section the excess Gibbs free energy 
models for binary solutions are presented. 

2.4.2. 1.1 Margules Equation 

The simplest expression for the excess Gibbs free energy for a binary solution is 
the two suffix Margules equation which is given by 

[g^ ! RT)= Ax,x^ (2.76) 

Where is an adjustable parameter which assume a constant value over a range 
of temperature and differ from system to system. Eq. (2.76) can be used to obtain the 
following equations relating the activity coefficients explicitly to liquid phase mole 
fractions. 

\nr,=Ax,^ iini) 

and ln/ 2 =^^i^ (2.78) 


2.4,2. 1.2 van Laar Equation 

The excess Gibbs free energy of a binary liquid solution can be represented by 
the two parameters van Laar Equation which is given by 



{g^IRT)-- 


Ax-^Xj 


X - 


hU/B) 


+ x. 


(2.79) 


The Eq. (2.79) can be manipulated to yield the following activity coefficient 
correlation. 



Iny - ^ 

[\-^(AIBlxJx,)f 

(2.80) 

and 

Iny - ^ 

[l+(5/4r,/x,)f 

(2.81) 


2.4.2. 1.3 Wilson Equation 


For mixtures in which the components differ from each other in molecular size 
and the interaction between the like and unlike molecules are different, Wilson 
proposed the following expression for the excess Gibbs free energy of a binary solution. 


(g^ /RT)= Xj ln(xj +F;2 xJ-X 2 ln(x2 +V2,x,) 


The Eq.(2.82) leads to the following activity coefficient correlations. 


and 


In^i =-ln(x, +F, 2 X 2 ) + x. 


ln^2 ~ ltl(x2 +1^21^1 



(2.82) 


(2.83) 

(2.84) 


2.4.2. 1.4 Non-Random Two Liquid (NRTL) Equation 


The NRTL equation proposed by Renon is applicable to partially miscible as 
well as completely miscible systems. The excess Gibbs free energy of a binary mixture 
is given by 




21^21 




'^12^12 


X | " I ” ^ 2^21 ^2 ^ 1^12 


/ Rt )= x , x , 


(2.85) 



Where 


^12 ={Si2 - 822)1 RT ( 2 . 86 ) 

^21 ={S2\~ Su)l RT (2.87) 

=exp(-a,3r,2) (2.88) 

<^21 =exp(-a,,T2j (2.89) 


The parameter cx .12 is introduced to take into account the non-randomness of the 
mixture. If Un^O, the mixture is completely random and equation Eq. (2.85) reduces to 
two suffix Margules equation. Usually cci 2 varies from 0.20 to 0.47 and in the absence 
of any information cell can be taken as 0.30. The parameter is introduced to take into 
account the energy of interaction between molecules of components i and j and it is 
analogous to Wilson’s parameter;!#,-. The activity coefficients can be derived from Eq. 
(2.85) and are given by 


and 


\nr, =^ 2 " 
In 72 


1 1 


1 ^12^12 

+ ^2^21 ^ 

(^2 +XjGi2) 

r,2^ 

' G., 1 

. 1 '^21 <^21 

J 

(xj H- X 2 G 21 ) 


(2.90) 

(2.91) 


2.4.3 Prediction of VLE data using UNIFAC method 


Fredenslund, Jones, and Pmsnitz first presented the UNIFAC GROUP 
contribution method. It was further developed by Fredenslund et. al., Gmehlig and 
Rasmussen, for use in practice. It is based on the concept that a liquid mixture may be 
considered as a solution of the stmctural units from which the molecules are formed 
rather then the solution of the molecules themselves (Concept given by Wilson and Deal 
in 1960). These structural units are called subgroups. To calculate the partial molar 
excess Gibbs free energies g,^ and thus to calculate the activity coefficients, the size 
parameter (volume and area parameters) and binary interaction parameter for each pair 
of subgroups are required. Size parameter will take care of shape and size of molecules 
and the interaction parameter gives the binary interaction between the two subgroups. 



The UNEFAC method for predicting liquid phase activity coefficients is 
essentially based on the UNIQUAC equation. In the UNIFAC method the activity 

coefficient consists of two parts - Iny,*^ (combinatorial) and In (residual). 


Iny,. =lny,.‘^ +Iny,.^ 


(2.92) 


The combinatorial part of the activity coefficient is given by 


and the residual part is given by 


tar,''=En"’(tar,-inr;") 


Where 


~ Qk 

Qm^m - 


V m J m 


= 




Surface area fraction of group m 


Wmn = exp 


^ u -un^ 

mn 


RT 


= expj 


^ a ^ 

mn 


V J- j 


)-('■! -l) 

^ = E v.'-'-'R, 

k 

?, =E''.“a 

k 

Qk = group volume parameter 
Rk = group area parameter 
Xm = mole fraction of group m in the mixture 
am„ = interaction parameter (in Kelvin ). ^nm 


(2.93) 


(2.94) 

(2.95) 

(2.96) 

(2.97) 

(2.98) 

(2.99) 

( 2 . 100 ) 


Umn = measure of interaction energy between groups m and n 



Fk = group residual coefficient and 

Fk^’^ = residual activity coefficient of group A: in a reference solution 
containing only molecules of type /. 

The parameters required for the calculation of activity coefficient using UNIFAC 
method are group area and volume parameters Qik, Qk) and the interaction parameter 
(um«) for each pair of subgroups. 



Chapter 3 

COMPUTER PROGRAM 


The software to calculate thennod 5 Tiamic properties contains, in all 77 programs, 
which are coupled by a program named MASTER program. Each program is responsible 
for a unique property prediction. Due to the excessive code length (Approximately 9000 
lines), it is not possible to provide source code for the software- instead simplified 
flowcharts are provided in this chapter. The programs for property estimation uses the 
methods described in Chapter 2. 

The problems of predicting the thermodynamic properties are divided in to several 
small parts to provide necessary flexibilily and to avoid the repeated coding for the same 
calculation. For example calculation of molar volume requires 

1 . The identification of EOS i.e. values of u and w, 

2. Evaluation of EOS constants a, b, A, B, S and a 

3 . Finally these u, w, a, b, A, B, S and a are used in the appropriate equation for the 
calculation of Molar Volume. 

The first two tasks are related with individual EOS, which are to be used for the 
calculation, and are common for almost all the programs. Hence, two different subroutines 
are developed (Get_UWJEOS for u and w values and Get_AB_EOS for calculating EOS 
parameters a, b, S and a). Similarly, all other properties calculations are divided into 
several small parts with the objective of reducing the code length and to provide flexibility. 

In this chapter flowcharts are provided for the individual property estimation 
methods followed by the flowcharts of all the subroutines, which are used in property 
estimation. 
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3.1 Calculation of Thermodynamic 
Properties ( Flowcharts ) 


Flowchart for the calculation of volume 



Fig. 3.1 Flowchart for volume calculation 








Flowchart for the calculation of Enthalpy and Entropy 

Departure 



Fig. 3.2 Flowchart for departure function calculation 
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Flowchart for the calculation ofFugacity 



Fig. 3.3 Flowchart forfugacity calculation 
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Flowchart for the calculation of Saturation Pressure 



Flowchart for Saturation pressure calculation 
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Flowchart for the calculation of Saturation Temperature 

C START ^ 



Fig. 3.5 Flowchart for Saturation temperature calculation 
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Flowchart for the calculation of Enthalpy and Entropy 



















A 


Call Entropy_Departure 

(Calculate Entropy Departure at 298. 15 K and ^ 29^5 • ^^ 2 ) 



Calculate 

h = h° + Ahi + Ah2 + Ahj + Ah4 
s = Asi + AS 2 + Ass + AS4 


Print h and s 


STOP 


Fig. 3.6 Flowchart for Enthalpy and Entropy calculation 













Flowchart for the calculation of Molar volume of Mixture 



Fig. 3.7 Flowchart for volume calculation of multi-component mixture 
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Flowchart for the calculation of Enthalpy and Entropy departure 

of Mixture 

( START ) 





Fig. 3.8 Flowchartfor Departure Junction calculation of multi-component mixture 
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Flowchart for the calculation of fugacity of Mixture 



Fig. 3,9 Flowchart for fugacity calculation of multi-component mixture 














Flowchart for the calculation of partial fugacity of a 
component in the Mixture 




Fig. 3.10 Flowchart for partial fugacity calculation of a component in a 

multi-component mixture 














Flowchart for the calculation of bubble point pressure 



Fig. 3.11 Flowchart for bubble point pressure calculation 














Flowchart for the calculation of bubble point temperature 



Fig. 3.12 Flowchart for bubble point temperature calculation 


















Flowchart for the calculation of Dew point temperature 



Fig 3.13 Flowchart for Dew point temperature calculation 
















Flowchart for the calculation of Dew point pressure 




Fig. 3.14 Flowchart for Dew point pressure calculation 















Flowchart for the calculation of T-x-y data 



Fig. 3.15 Flowchart for VLB calculation 

















Subroutine Get_UW_EOS 



Fig. 3.16 Subroutine for u and w 
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Subroutine Get AB EOS 



Fig. 3.17 Subroutine for EOS parameter 
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Subroutine Volume 



Fig. 3.18 Subroutine for volume calculation 
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Subroutine Enthalpy Departure 



Fig. 3.19 Subroutine for enthalpy departure calculation 


Subroutine Entropy Departure 



Fig. 3.20 Subroutine for Entropy departure calculation 











Subroutine Fugadiy 



Fig. 3,21 Subroutine for fugacity calculation 


Subroutine Watson 
















Fig. 3.22 Subroutine for Watson correlation 


Subroutine NTPNTHJDL 



Fig. 3.23 Subroutine for enthalpy change for ideal gas 
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Subroutine vdWMixing 



Fig. 3.24 Subroutine for vdWMixing rule 

Subroutine PartialFugacity 



Fig. 3.25 Subroutine for Partial fugacity 
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Subroutine Margules 



Fig. 3.26 Subroutine for Margules Equation 


Subroutine vanLaar 



Fig. 3.27 Subroutine for van Laar Equation 
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Subroutine Wilson 



Fig. 3.28 Subroutine for Wilson Equation 

Subroutine NRTL 



Fig. 3.29 Subroutine for NRTL Equation 
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Subroutine UNIFAC 



Fig. 3.30 Subroutine for UNIFAC Equation 


Lt- 
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Subroutine Sat T 



Fig. 3.31 Subroutine satJT 
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Subroutine Sat_P 



Fig. 3.32 Subroutine sat P 













Chapter 4 

SAMPLE INPUT / OUTPUT 


In this chapter, sample input/output is provided for each of the thermodynamic 
properties that can be predicted by the software. As the input windows are same for 
almost all the properties, to avoid the repetition the input windows are shown for 
volume calculation only. For the rest of the properties only the output window is 
provided which contains input as well as output information. 

4.1 Sample Input/Output For Volume Calculation 

Here is a sample input/output for the calculation of pure component volume 
using cubic equation of state. The first window that appears is shown below. 



Volume Calculation using Cubic EOS 


Volume Calculation 



Set Units for T, P, & V 
Run Volume Calculator 
Quit 





When user selects Run volume calculator and presses the ENTER Key the 
following input window will appear in which the user enters the values shown in italics. 




VOLUME CALCULATION USING EOS 


The CARDAN’S Method will be used To Solve Cubic EOS. This method allows 
direct solution of any cubic Equation Of State. The following four cubic EOS can be 
selected. 


1 . van der Walls (VDW) 2. Redlich Kwong (RK) 

3. Soave Redlich Kwong (SRK) and 4. Peng Robinson (TR) EOS. 


Enter the Temperature of the Fluid in K 

; 427.98 

Enter the Acentric Factor 'w' 

; 0.298 

Type of Equation to use? VDW/RKTPR/SRK 

: PR 



When user presses FI Key after entering the above data, new Window will 
appear asking for the name of substance for which calculation is to be done. 


Enter The Name of the Substance: N-OCTANE 


This completes the required input for the program and the output window will 
appear as shown below. 







Volume Calculations using PENG ROBINSON equation of state 


Input Data: 

Compound Name: N-OCTANE 

Critical Temperature = 569.8 K Critical Pressure = 24.970 BAR 

Temperature = 428.0 K Pressure = 2.150 BAR 

Acentric Factor = 0.3980 


All the Calculations were done in the following units 

Temperature in Kelvin Pressure in Pascal Volume in m3 /mol. 

Input Data ( After Conversion) 


Critical Temperature = 569.8 Kelvin Critical Pressure = 2497000 Pascal 
Temperature = 428.0 Kelvin Pressure = 215000 Pascal 

Acentric Factor = 0.398000 


A =0.088494 B =0.008919 

alpha = -0.99 1081 beta = 0.070418 gamma = -0.000709 

p =-0.256996 q =-0.049555 D =-0.000015 


The value of theta is 0. 1 53 645 
Vg = 1.51 4e-02 cubic meter/mol 
VI = 2.003e-04 cubic meter/mol 

RESULTS: 


Vg( per mole) = 1 5 . 143 LITERS 
Vl( per mole) = 0.200 LITERS 



The Volume calculator can give P-v as well as T-v diagram for which the user will 
provide temperature and pressure range. Following is a P-v diagram for T= 450.0 
Kelvin and for the pressure range of I.O bar to 6.0 bar. 





Uolume in riB 


Fig. 4.1. Pressure vs. Molar volume diagram 

4.2 Sample Output For Enthalpy And Entropy 
Departure Calculation 

Following is the output window, obtained by using the departure calculator to 
calculate enthalpy and entropy departure for N-OCTANE at 428. 0 Kelvin and 2. 15 bar. 



Departure Calculations for PENG ROBINSON equation of state 

Input Data; 



Compound Name; N-OCTANE 
Critical Temperature = 569. 8 K 

Critical Pressure = 

24.970 BAR 

Temperature = 428.0 K 

Acentric Factor = 0.3980 

Pressure = 

2.150 BAR 



Input Data ( After Conversion) 


Critical Temperature = 569.8 Kelvin Critical Pressure = 
Temperature = 428.0 Kelvin Pressure 

Acentric Factor = 0.398000 


2497000 Pascal 
215000 Pascal 







Ail the Calculations were done in the following units 
T emperature in Kelvin Pressure in Pascal 


Volume in m3/moL 


A =0.088325 B =0.008912 

Alpha =-0.991088 beta = 0.070262 gamma = -0.000707 

p =-0.257157 q =-0.049606 D =-0.000015 

The value of theta is 0.153058 
Vg = 1 .5 15e-02 cubic meter/mol 
VI = 2.003e-04 cubic meter/mol 

RESULTS: 


dh(ENTHALPY DEPARTURE) = -0.8896 kJ/mol 
ds(ENTROPY DEPARTURE) = -1.3973 J/mol.K 



The Departure calculator can give depctrture function vs'. pressure diagram as 
well as departure function vs. temperature diagram for which the user will provide 
temperature and pressure ranges. The following enthalpy departure vs'. pressure 
diagram and entropy departure vs pressure diagram are obtained by using departure 
calculator for N-OCTANE at T= 427.0 Kelvin and in the pressure range of 1.0 bar to 
IQ.Obar. 



0.8 1.8 2.8 3.8 4.8 5.8 6.B 7.8 

Pre-ss-ure in 127.92857 kpa 


Fig. 4.2 Enthalpy departure vs. pressure diagram 






Pressiure In 12S.G428G kpa 


Fig. 4.3 Entropy departure vs. pressure diagram 

4.3 Sample Output For Fugacity Calculation 

The output window obtained by using the fugacity calculator for N-OCTANE at 427.98 
Kelvin and 2. 15 bar is shown below. 


mm 




Fugacity Calculations using PENG ROBINSON equation of state 

Input Data; 



Compound Name: N-OCTANE 

Critical Temperature = 569.8 K 

Temperature = 428.0 K 

Acentric F actor = 0.3980 

Critical Pressure = 
Pressure = 

24.970 BAR 

2.150 BAR 

All the Calculations were done in the following units 

Temperature in Kelvin Pressure in Pascal 

Volume in m3 /mol. 

Input Data ( After Conversion) 



Critical Temperature = 569.8 Kelvin 

Critical Pressure = 

2497000 Pascal 

Temperature = 428.0 Kelvin 

Pressure = 

215000 Pascal 

Acentric Factor = 0.398000 










A = 0.088494 
alpha = -0.991081 
p = -0.256996 


B =0.008919 
beta = 0.070418 
q = -0.049555 


gamma = -0.000709 
D =-0.000015 


The value of theta is 0.153058 
Vg = 1.515 e-02 cubic meter/mol 
VI = 2. 003 e-04 cubic meter/mol 
Fugacity coefficient of liquid = 
Fugacity coefficient of vapor = 


0.909 

0.921 


RESULTS; 

FUGACITY OF LIQUID 
FUGACITY OF VAPOR 


= 1.9552 BAR 

= 1.9808 BAR 



The fugacity calculator can give fugacity coefficient V5. pressure as well as 
fugacity coefficient vs temperature diagram for which the user will provide temperature 
and pressure range. Following is a fugacity coefficient vs pressure diagram for T= 450. 0 
Kelvin and for the pressure range of 1. 0 bar to 6. 0 bar. 



Fig. 4.4 Fugacity coefficient of liquid vs. pressure diagram 







2.8 3.8 4.8 5.3 &.8 7.B 8 8 3.8 


Pressure in 35.35714 kPa 


Fig. 4.5 Fugacity coefficient of vapor vs. pressure diagram 

4.4 Sample Output For Saturation Temperature 
Calculation 

Following is the output window, obtained by using the saturation temperature 
calculator to calculate Saturation temperature ^oxN-OCTANE at 2.15 bar. 



Sat. Temperature Calculations using PENG ROBINSON equation of state 
Input Data: 

Compound Name; N-OCTANE 

Critical Temperature = 569.8 K Critical Pressure = 24.970 BAR 

Pressure = 2. 150 BAR 

Acentric Factor = 0.3980 


All the Calculations were done in the following units 

Temperature in Kelvin Pressure in Pascal Volume in m3/mol 

Input Data ( After Conversion) 

Critical Temperature = 569.8 Kelvin Critical Pressure = 2497000 Pascal 
Pressure = 215000 Pascal 

Acentric Factor = 0.398000 





The saturation temperature calculator can give the following saturation 
temperature vs. pressure diagram 



Fig. 4.6 Saturation temperature V5. pressure diagram 


4.5 Sample Output For Saturation Pressure 
Calculation 

Following is the output window obtained by using saturation pressure calculator 
foxN-OCTANE at 450.0 Kelvin. 





L 


Sat. Pressure Calculations using PENG ROBINSON equation of state 


Input Data: 

Compound Name: N-OCTANE 
Critical Temperature = 569.8 K 
Temperature = 450.0 K 

Acentric Factor = 0.3980 


Critical Pressure = 24.970 BAR 


All the Calculations were done in the following units 

T emperature in Kelvin Pressure in Pascal Volume in m3/mol 

Input Data ( After Conversion) 

Critical Temperature = 569.8 Kelvin Critical Pressure = 2497000 Pascal 
Temperature = 450.0 Kelvin 

Acentric Factor = 0.398000 


Iteration no:- 1 
Guess Pressure = 
Iteration no:- 2 
Guess Pressure = 


347.1 kPa Phi L= 0.889 


346.4 kPa 


Phi L= 0.891 


Phi_v= 0.891 
Phi v= 0.891 


RESULT: 

Saturation Pressure = 


3.5 BAR 







4.6 Sample Output For Vapor Enthalpy 
Calculation 

Following is the output window, obtained by using the Enthalpy calculator to 
calculate vapor enthalpy for N-OCTANE at 428. 0 Kelvin and 2. 15 bar. 



Vapor Enthalpy Calculations using PENG ROBINSON equation of state 




Hv at NBP 

dhl (Hv at 298.15 K (1.976 kPa ) ) 

dh2 (Enthalpy departure) At 298. 15 K And 1.976 kPa 

dhS (Enthalpy departure) At 427.980 K And 215.0 kPa 


= 34.864kJ/niol dh 
= 41.603 kJ/mol 
= 0.013 kJ/mol 
= -0.890 kJ/mol 


dh4 (Enthalpy change due to) Change in the ideal gas state at 298.15 K and 
1.976 kPa to ideal gas state at 427.980 K and 215.000 kPa = 21.954 kJ/mol 

RESULT: 

ENTHALPY OF VAPOR = 65. 128 kJ/mol 





The Software can give the following Enthalpy vs. Pressure rliagrqr>i for N-OCTANE 
stT= 450 Kelvin and in the Pressure range of 1. 0 bar to 5.0 bar. 
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Fig. 4.8 Enthalpy vx Pressure Diagram 
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4.7 Sample Output For Vapor Entropy 
Calculation 

Following is the output window, obtained by using the Entropy calculator to 
calculate vapor enthalpy ^or N-OCTANE at 428.0 Kelvin and 2.15 bar. 



Vapor Entropy Calculations using PENG ROBINSON equation of state 

Input Data: 


Compound Name: N-OCTANE 


Critical Temperature = 569.8 K 
Temperature ' = 428.0 K 

Acentric Factor = 0.3980 

Normal boiling point (NBP) = 398.8 K 

Critical Pressure = 24.970 BAR 

Pressure = 2.15 BAR 


All the Calculation were done in the following units 
Temperature in Kelvin Pressure in Pascal 


Volume in m3 /mol 







Input Data ( After Conversion) 



Critical Temperature 

= 569.8 Kelvin 

Critical Pressure 

= 2497000 Pascal 

Temperature 

= 428.0 Kelvin 

Pressure 

= 215000 Pascal 

Acentric Factor 

= 0.398000 

NBP 

= 398.8 Kelvin 


Hv at NBP = 34. 864kJ/mol 

dsOEntropy change of vaporization)AT 298. 1 5 K( 1 .976 kPa) = 139.537 J/mol. 
ds(Entropy departure) At 298. 150 K And 1.976 kPa = 0.027 J/mol. 

ds(Entropy departure) At 428.000 K And 215.000 kPa = -1.397 J/mol. 

ds(Entropy change due to)Change in the ideal gas 
state at 298.15 K and 1.976 kPa to ideal gas 

state at 428.000 K and 215.000 kPa = 21.843 J/mol. 

RESULT: 

ENTROPY OF VAPOR = 159.96 J/mol.K 








4.8 Sample Output For Molar Volume 
Calculation of Mixture 

Following is the output window, obtained by using the Software to calculate 
Molar Volume of N-Octane-Toluene mixture at 600.0 Kelvin and 16.0 bar. Mole 
fraction ofN-Octane = 0.4 



Calculation of Molar Volume of a Multi-component Mixture 


Input Data: 


Component Name 

Tc(K) 

Pc (BAR) 

W 

Y 

N-OCTANE [1] 

568.8 

24.97 

0.3980 

0.4000 

TOLUENE [2] 

591.8 

41.00 

0.2630 

0.6000 


Mixture Temperature and Pressure;- 

T(K)= 600.0 P(BAR)= 16.000 

EOS used For Calculation is PR 

For the calculations Unit of T is taken as Kelvin and Unit of P is Pascal 


OUTPUT 

Component No 

1 

2 

Mixture 


a (Pa/moU-2) 

3.8991 

2.6715 

3.1348 


b (m'^3/mol) 

1.48e-04 

9.34e-05 

1.15e-04 


For the Mixture parameter (a & b) calculation van der Waals mixing rule had 
been used. 


A =0.201559 
alpha = -0.963075 
p =-0.185552 


B = 0.036925 
beta = 0.123620 
q =-0.032512 


gamma = -0.006029 
D = 0.000028 


Vg 


2.561e-03 cubic meter/mol 



RESULT 

Molar Volume of Vapor Mixture = 2.561e-03 mS 



4.9 Sample Output For Enthalpy And Entropy 
Departure Calculation of Mixture 

Following is the output window, obtained by using the Software to calculate 
enthalpy and entropy departure of Ethanol-Methanol mixture at 450.0 Kelvin and 10.0 
bar. Mole jfraction oi Methanol -0.4 



Calculation of Departure functions of a Multi-component Mixture 


Input Data: 


Component Name 

Tc(K) 

Pc (BAR) 

W 

Y 

ETHANOL [1] 

513.9 

61.40 

0.6640 

0.6000 

METHANOL [2] 

512.6 

80.90 

0.5560 

0.4000 


Mixture Temperature and Pressure; - 


T(K)= 450.0 P(BAR)= 10.000 

EOS used For Calculation is PR 


For the calculations Unit of T is taken as Kelvin and Unit of P is Pascal 


OUTPUT 

Component No a (Pa/moU-2) b (m^3/mol) 


1 

2 

Mixture 


1.5921 5.41e-05 

1.1806 4.10e-05 

1.4201 4.89e-05 






For the Mixture parameter (a & b) calculation van der Waals mixing rule had 
been used 


Vapor compressibility factor ( mixture) = 0.9057 

RESULT 


Enthalpy Departure of Mixture = -1.210 kJ/mol 

Entropy Departure of Mixture = - 1 . 927 J/mol.K 



4.10 Sample Output For Fugacity Calculation of 
Mixture 

Following is the output window, obtained by using the Software to calculate 
Fugacity of Benzene-Toluene mixture at 600.0 Kelvin and 16.0 bar. Mole fraction of 
Toluene = 0.5. 



Calculation of fugacityof a Multi-component Miiture 


Input Data: 


Component Name 

Tc(K) 

Pc (BAR) 

W 

Y 

BENZENE [1] 

562.2 

48.90 

0.2120 

0.5000 

TOLUENE [2] 

591.8 

41.00 

0.2630 

0.5000 


Mixture Temperature and Pressure;- 

T(K)= 600.0 P(BAR)= 16.000 


EOS used For Calculation is PR 

For the calculations Unit of T is taken as Kelvin and Unit of P is Pascal 






OUTPUT 

Component No a (Pa/mol'^-2) 

b (m^3/mol) 

1 1.9508 

2 2.6715 

Mixture 2.2970 

7.44e-05 

9.34e-05 

8.39e-05 

For the Mixture parameter (a & b) calculation van der Waals mixing rule 
had been used. 

A= 0.148 B = 2.690e-02 Z(Liquid) = 0.8722 

Z(Vapor)= 0.8722 

Vapor : PHI = 0.8833 Fugacity = 1413.265 kPa 

RESULT 


Fugacity of vapor Mixture = 14.133 BAR 


1 ^— 

■raBB 


4.11 Sample Output For Partial Fugacity 
Calculation of Mixture 

Following is the output window, obtained by using the Software to calculate 
Partial Fugacity of Benzene-Toluene mixture at 600.0 Kelvin and 16.0 bar. Mole 
fraction of Benzene = 0.5 


Calculation of Partial Fugacity a Multi-component Mixture 


Input Data: 


Component Name 

Tc(K) 

Pc (BAR) 

W 

Y 

BENZENE [1] 

562.2 

48.90 

0.2120 

0.5000 

TOLUENE [2] 

591.8 

41.00 

0.-2630 

0.5000 



Mixture Temperature and Pressure:- 








T(K)= 600.0 P(BAR)= 16.000 

EOS used For Calculation is PR 

For the calculations Unit of T is taken as Kelvin and Unit of P is Pascal 


OUTPUT 

Component No a (?dJmoV^-2) b (m'^3/mol) 

1 1.9508 7.44e-05 

2 2.6715 9.34e-05 

Mixture 2.3696 8.58e-05 

For the Mixture parameter (a & b) calculation van der Waals mixing rule 
had been used 


For Vapor : A= 0.1524 B = 2.75e-02 Z= 0.8674 

phi[l]= 0.9031 f[l]= 577.983 kPa 

phi[2]= 0.8640 f[2]= 829.440 kPa 

RESULT 


Partial Fugacity of Component 1 (Vapor) = 5.780 BAR 
Partial Fugacity of Component 2 (Vapor) = 8.294 BAR 



4.12 Sample Output For Bubble Point Pressure 
Calculation 

Following is the output window, obtained by using the Software to calculate 
Bubble Point Pressure of Ethane-Propane mixture at 200.0 Kelvin. Mole fi'action of 
Ethane = 0.7 


""1 



Calculation of Bubble Pressure of a Multi-component Mixture 


Input Data: 

Component Name Tc(K) Pc (BAR) W 


ETHANE [1] 
PROPANE [2] 


305.4 

369.8 


48.80 

42.50 


0.0990 0.7000 

0.1530 0.3000 




Mixture T emperature: - T (K) = 200. 0 

EOS used For Calculation is PR 

For the calculations Unit of T is taken as Kelvin and Unit of P is Pascal 
OUTPUT: ~~ “■ 

Calculation of bubble pressure : 

Saturation Pressure of indeveduaJ Component At T = 200 0 K 


SATP[I] 
SAT P [2] 


218 kPa 
21 kPa 


Itr No P guess 

Y[l] 

Y[2] 

K[l] 

K[2] 

I 

79762.046 

0.79430 

0.18618 

2.6477 

0.2660 

1 

79762.046 

0.79430 

0.18618 

2.6477 

0.2660 

2 

79762.046 

0.79434 

0.18614 

2.6478 

0.2659 

3 

79762.046 

0.79434 

0.18614 

2.6478 

0.2659 

5 

78204.982 

0.80971 

0.18980 

2.6990 

0.2711 

6 

78204.982 

0.80977 

0.18975 

2.6992 

0.2711 

7 

78204.982 

0.80977 

0.18975 

2.6992 

0.2711 

9 

78167.445 

0.81009 

0.18989 

2.7003 

0.2713 

10 

78167.445 

0.81014 

0.18984 

2.7005 

0.2712 

11 

78167.445 

0.81014 

0.18984 

2.7005 

0.2712 

13 

78166.522 

0.81010 

0.18989 

2.7003 

0.2713 

14 

78166.522 

0.81015 

0.18985 

2.7005 

0.2712 

15 

78166.522 

0.81015 

0.18985 

2,7005 

0.2712 


RESULT: 


Bubble Pressure of The Mixture = 
Composition of the Mixture; 
Component Name Y 

ETHANE 0.8102 
PROPANE 0.1898 


= 0.782 BAR 




4. 1 3 Sample Output For Bubble Point 
Temperature Calculation 

Following is the output window, obtained by using the Software to calculate 
Bubble Point Temperature of Ethane-Propane mixture at 15.0 bar. Mole fraction of 

Ethane ~ 0.3 



Calculation of Bubble Temperature of a Multi-component Mixture 


Input Data: 


Component Name 

Tc(K) 

Pc (BAR) 

W 

Y 

ETHANE [I] 

305.4 

48.80 

0.0990 

0.3000 

PROPANE [2] 

369.8 

42.50 

0.1530 

0.7000 

Mixture Pressure 

P (Bar) = 

15.0 




EOS used For Calculation is PR 

For the calculations Unit of T is taken as Kelvin and Unit of P is Pascal 


OUTPUT 

Calculation Of Bubble Temperature: 
P = 1500.0 kPa 


SATT[1] = 
SAT T [2] = 

Itr No T guess 

255 Kelvin 

3 17 Kelvin 
Y[l] 

Y[2] 

K[l] 

K[2] 

1 

282.188 

0.67414 

0.38757 

1.6854 

0.6460 

2 

283.188 

0.64689 

0.31248 

1.6172 

0.5208 

3 

284.188 

0.66530 

0.31290 

1.6632 

0.5215 

4 

285.188 

0.67442 

0.32092 

1.6861 

0.5349 

5 

284.188 

0.68334 

0.32908 

1.7083 

0.5485 

6 

284.288 

0.67246 

0.32303 

1.6811 

0.5384 

7 

284.388 

0.67448 

0.32275 

1.6862 

0.5379 

8 

284.488 

0.67538 

0.32354 

1.6884 

0.5392 

9 

284.388 

0.67626 

0.32435 

1.6906 

0.5406 



10 

11 

12 

13 

14 

15 

16 

17 

18 

19 

20 
21 
22 

23 

24 

25 


284,398 

284.408 

284.418 

284.428 

284.438 

284.448 

284.458 

284.448 

284.449 

284.450 

284.451 

284.452 

284.453 
284.452 
284.452 
284.452 


RESULT: 
Tbubble (In K ) 


0.67518 

0.67537 

0.67546 

0.67555 

0.67564 

0.67573 

0.67581 

0.67590 

0.67579 

0.67581 

0.67582 

0.67583 

0.67584 

0.67585 

0.67584 

0.67584 


0.32374 

0.32371 

0.32379 

0.32387 

0.32396 

0.32404 

0.32412 

0.32420 

0.32414 

0.32413 

0.32414 

0.32415 

0.32416 

0.32417 

0.32416 

0.32416 


1.6879 

1.6884 

1.6887 

1.6889 

1.6891 

1.6893 

1.6895 

1.6898 

1.6895 

1.6895 

1.6896 
1.6896 
1.6896 
1.6896 
1.6896 
1.6896 


284.45 Y[l]= 0.67584 Y[2]= 0.32416 


0.5396 

0.5395 

0.5397 

0.5398 

0.5399 

0.5401 

0.5402 

0.5403 

0.5402 

0.5402 

0.5402 

0.5403 

0.5403 

0.5403 

0.5403 

0.5403 



J: 

j 

4. 1 4 Sample Output For Dew Point Pressure 
Calculation 

Following is the output window, obtained by using the Software to calculate 
Dew Point Pressure of Ethane-Propane mixture at 200.0 Kelvin. Mole fraction of 
Ethane = 0.3. 



Calculation of Dew Pressure of a Multi-component Mixture 


Input Data: 


Component Name 

Tc(K) 

Pc (BAR) 

W 

Y 

ETHANE [1] 

305.4 

48.80 

0.0990 

0.3000 

PROPANE [2] 

369.8 

42.50 

0.1530 

0.7000 




Mixture Temperature:- T (K) = 200,0 
EOS used For Calculation is PR 

For the calculations Unit of T is taken as Kelvin and Unit of P is Pascal 


OUTPn': 


Calculation Of Dew Pressure : 

Saturation Pressure of individual Component At 
T = 200.0 K 


SAT P [13 = 217761 kPa 




SAT P [2] = 20620 kPa 




Itr No 

Pguess(Pa) 

X[13 

X[2] 

K[l] 

K[2] 

1 

79762.046 

0.10233 

2.44851 

2.9318 

0.2859 

2 

79762.046 

0.10311 

2.46667 

2.9095 

0.2838 

3 

79762.046 

0.10305 

2.46528 

2.9112 

0.2839 

4 

79762.046 

0.10306 

2.46538 

2.9110 

0.2839 

5 

79762.046 

0.10306 

2.46538 

2.9110 

0.2839 

6 

79762.046 

0.10306 

2.46538 

2.9110 

0.2839 

7 

79762.046 

0.10306 

2.46538 

2.9110 

0.2839 

8 

79762.046 

0.10306 

2.46538 

2.9110 

0.2839 

10 

31054.768 

0.04414 

1.04806 

6.7959 

0.6679 

11 

31054.768 

0.04416 

1.04821 

6.7928 

0.6678 

12 

31054.768 

0.04416 

1.04821 

6.7929 

0.6678 

13 

31054.768 

0.04416 

1.04821 

6.7928 

0.6678 

14 

31054.768 

0.04416 

1.04821 

6.7928 

0.6678 

16 

28428.785 

0.04077 

0.96203 

7.3582 

0.7276 

17 

28428.785 

0.04078 

0.96203 

7.3570 

0.7276 

18 

28428.785 

0.04078 

0.96203 

7.3570 

0.7276 

20 

28349.336 

0.04067 

0.95937 

7.3769 

0.7296 

21 

28349.336 

0.04067 

0.95937 

7.3757 

0.7296 

22 

28348.177 

0.04067 

0.95933 

7.3771 

0.7297 

25 

28348.177 

0.04067 

0.95933 

7.3760 

0.7297 

26 

28348.177 

0.04067 

0.95933 

7.3760 

0.7297 



RESULT: 

Dew Point Pressure of The Mixture = 0.283 BAR 

Composition of the Mixture: 

Component Name X 
ETHANE 0.04067 
PROPANE 0.95933 



4.15 Sample Output For Dew Point 
Temperature Calculation 

Following is the output window, obtained by using the Software to calculate 
Dew Point Temperature of Ethane-Propane mixture at 13.0 bar. Mole fraction of Ethane 



Calculation of Dew Temperature of a Multi-component Mixture 


Input Bata; 
Component Name 

Tc(K) 

Pc (BAR) 

W Y 

ETHANE [1] 
PROPANE [2] 

305.4 

369.8 

48.80 

42.50 

0.0990 0.3000 

0.1530 0.7000 

Mixture Pressure 

P (Bar) = 

13.0 


EOS used For Calculation is PR 

For the calculations Unit of T is taken as Kelvin 

and Unit of P is Pascal 


OUTPUT: 

Calculation Of Dew Temperature: 

Saturation Temperature of individual Component At 


p 


1300.0 kPa 


Itr No T guess 

Y[l] 

Y[2] 

K[l] 

K[2] 


! 

290.000 

0.14680 

1.03868 

2,0436 

0.6739 


2 

300.000 

0.12894 

0.84933 

2.3267 

0.8242 


3 

290.000 

0.14680 

1.03868 

2.0436 

0.6739 


4 

291.000 

0.14478 

1.01718 

2.0722 

0.6882 


5 

292.000 

0.14281 

0,99630 

2.1007 

0.7026 


6 

293.000 

0.14089 

0.97602 

2.1293 

0.7172 


7 

294.000 

0.13903 

0.95633 

2.1578 

0.7320 


8 

295.000 

0.13723 

0.93720 

2.1862 

0.7469 


9 

296.000 

0.13547 

0.91861 

2.2145 

0.7620 


10 

297.000 

0.13376 

0.90054 

2.2428 

0.7773 

1 

11 

298.000 

0.13211 

0.88299 

2.2709 

0.7928 

■1 

12 

299.000 

0.13050 

0.86592 

2.2989 

0.8084 

4 

13 

298.000 

0.13211 

0.88299 

2.2709 

0.7928 


14 

298.100 

0.13194 

0.88126 

2.2737 

0.7943 

1 

15 

298.200 

0.13178 

0.87954 

2.2765 

0.7959 

.1 

16 

298.300 

0.13162 

0.87782 

2.2793 

0.7974 


17 

298.400 

0.13146 

0.87610 

2.2821 

0.7990 


18 

298.500 

0.13130 

0.87439 

2.2849 

0.8006 

1 

: 

19 

298.600 , 

0.13114 

0.87269 

2.2877 

0.8021 


20 

298.700 

0.13098 

0.87099 

2.2905 

0.8037 


21 

298.800 

0.13082 

0.86930 

2.2933 

0.8052 


22 

298.900 

0.13066 

0.86761 

2.2961 

0.8068 

!' 

23 

298.800 

0.13082 

0.86930 

2.2933 

0.8052 


24 

298.810 

0.13080 

0.86913 

2.2936 

0.8054 


25 

298.800 

0.13082 

0.86930 

2.2933 

0.8052 


26 

298.801 

0.13082 

0.86928 

2.2933 

0.8053 


27 

298.802 

0.13081 

0.86926 

2.2933 

0.8053 


28 

298.803 

0.13081 

0.86925 

2.2934 

0.8053 


29 

298.804 

0.13081 

0.86923 

2.2934 

0.8053 


30 

298.805 

0.13081 

0.86921 

2.2934 

0.8053 


31 

298.806 

0.13081 

0.86919 

2.2934 

0.8053 


32 

298.807 

0.13081 

0.86918 

2.2935 

0.8054 


33 

298.806 

0.13081 

0. 86919 

2.2934 

0.8053 




34 298.806 


0.13081 


0.86919 2.2935 


0.8053 


RESULT: 

Tdew(InK)= 298.81 X[l] = 0.13081 X[2]= 0.86919 



4. 1 6 Sample Output For VLB data Calculation 
using Excess Gibbs Free Energy Models 

Following is the output window, obtained by using the Software to calculate 
VLB (T-x-y) data for Ethanol-Benzene mixture at 1.0 bar. 



Calculation of VLE (T-x-y) data for binary mixture 


Input Data: 
Component Name 

Tc(K) 

Pc (BAR) 

W 

ETHANOL [1] 

513.9 

61.40 

0.6640 

BENZENE [2] 

562.2 

48.90 

0.2120 

Mixture Pressure:- 

P(BAR) = 

1.0 



EOS used For Calculation is PR 
Acti’vity coefficient model : van Laar 

For the calculations Unit of T is taken as Kelvin and Unit of P is Pascal 


OUTPUT: 

Pressure ( At which VLE data is required ) = 100.000 kPa 





Component Saturation Tempereture 

1 351.49 

2 352.95 

T-X-Y data 


No Of 
iteration 

T(K) 

X[l] 

Y[l] 

Gamma[l] 

Gamma[2] 

105 

352.8 

0.00000 

0.00000 

6.4045 

1.0000 

97 

343.8 

0.10000 

0.30616 

4.1788 

1.0224 

93 

341.5 

0.20000 

0.39028 

2.9315 

1.0881 

98 

340.8 

0.30000 

0.42526 

2.1919 

1.1985 

94 

340.6 

0.40000 

0.44473 

1.7334 

1.3594 

89 

340.6 

0.50000 

0.46183 

1.4403 

1.5813 

97 

340.8 

0.60000 

0.48504 

1.2502 

1.8794 

97 

341.4 

0.70000 

0.52342 

1.1277 

2.2748 

89 

342.8 

0.80000 

0.59162 

1.0525 

2.7965 

87 

345.7 

0.90000 

0.72123 

1.0123 

3.4832 

97 

351.4 

1.00000 

1.00000 

1.0000 

4.3864 



4.17 Sample Output For VLE data Calculation 
using UNIFAC Models 

Following is the output window, obtained by using the Software to calculate 
VLE (T-x-y) data fox Methanol-Water mixture at 1.0 bar. 



Calculation of VLE (T-x-y) data for binary mixture 

VLE CALCULATION FOR METHANOL [1] -WATER [2] SYSTEM 
Name of the subgroups present in the binary mixture and number of times 



They appear in the individual component. 



Name Of Sub-Groups 

METHANOL 

WATER 

CH30H 



1 


0 

H20 



0 


1 

Pressure = 

LOBAR 





3|c5fc5ie5ie5ic5le%3jc5je^3|c5|c:ic:ic5(c3jc:(c:)|c%5(c5lc 

VLE-DATA ********^************* 

T(0C) 

XI 

Y1 


acl 

ac2 

♦ >|c*******j|c* + ****^£:)!**:)cj|c:(ej|:i|c*:)cHc****))c**^***,|c:(c*j|e**!|c*:)c*J(c»4:* 

99.11 

0.000 

0.000 


1.8431 

1.0000 

88.73 

0.100 

0.383 


1.5348 

1.0094 

83.02 

0.200 

0.551 


1.3439 

1.0331 

79.11 

0.300 

0.653 


1.2210 

1.0666 

76.06 

0.400 

0.728 


1.1396 

1.1069 

73.49 

0.500 

0.788 


1.0849 

1.1524 

71.22 

0.600 

0.840 


1.0483 

1.2017 

69.16 

0.700 

0.885 


1.0245 

1.2542 

67.24 

0.800 

0.926 


1.0099 

1.3092 

65.45 

0.900 

0.965 


1.0023 

1.3664 

63.74 

1.000 

1.000 


1.0000 

1.4253 

3^ 3j|C 3ljc 3^ jjc )?||C ]3||C Sjc 3^ IjjC iSjC 3fjc ijc 3|c 3^ 3jc p|C Sjc 5|c ijJC 3jC 3||C 5jc sjC sjc Sjc 3jc sfC SjjC Sfc sfc jjc jjc 3^ 3^ 3|c 3|c sjc 3^ 3jc sjf 3|C 








The Software can give T-x-y and x-y diagrams. The T-x-y and x-y diagrams 
obtained for Methanol-Water system at P = / bar is given in Figs. 4.10 and 4.11, 
respectively. 
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Fig. 4.10 T-x-y diagram for Methanol [1]-Water [2] system at 1 bar. 
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Fig. 4.11 x-y diagram for Methanol [lJ-Water[2J system at 1 bar. 
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Chapter 5 


CONCLUSIONS AND 
SUGGESTIONS FOR FURTHER 

WORK 


This software is especially designed for the estimation of thermodynamic 
properties. It is user friendly and requires no prior knowledge of any computer 
language. As the software has been developed on MS-DOS operating system, it can be 
used only with computers, which have MS-DOS platform. 

The software is a step in the direction of the development of an extensive 
software for chemical engineering applications. There is a need to include some more 
properties like viscosity, thermal conductivity, difiusion coefficients etc. For more 
accurate predictions newer equations of states can be used with some of the recently 
developed miring rules. As for the phase equilibrium calculations liquid-liquid 
equilibrium, solid-liquid equilibrium etc. can be added. 
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